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ABSTRACT 

We show that radiation dominated accretion discs are likely to suffer from a "pho- 
ton bubble" instability similar to that described by Arons in the context of accretion 
onto neutron star polar caps. The instability requires a magnetic field for its existence. 
In an asymptotic regime appropriate to accretion discs, we find that the overstable 
modes obey the remarkably simple dispersion relation 

uj 2 = -igkF(B,k). 

Here g is the vertical gravitational acceleration, B the magnetic field, and F is a geo- 
metric factor of order unity that depends on the relative orientation of the magnetic 
field and the wavevector. In the nonlinear outcome it seems likely that the instabil- 
ity will enhance vertical energy transport and thereby change the structure of the 
innermost parts of relativistic accretion discs. 
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1 INTRODUCTION 

Compact objects such as black holes and neutron stars are 
among the most interesting objects in astrophysics because 
of their exotic strong-field gravitational physics and because 
they are likely at the center of some of the most luminous 
and energetic objects in the universe. Since it is mainly their 
accretion flows, and not the compact objects themselves, 
that are readily observed, the dynamics and radiative prop- 
erties of the accretion flows are the focus of much theoretical 
attention. 

One of the most important models for accretion flows 
is the thin disc (Shakura & Sunyaev 1973, Lynden-Bell & 
Pringle 1974). This model predicts that for rapidly accret- 
ing objects the inner parts of the disc are radiation pressure 
dominated. It was realized soon after the model was first 
proposed, however, that radiation pressure dominated discs 
are unstable, both viscously (Lightman & Eardley 1974) 
and thermally (Shakura & Sunyaev 1975). This suggests 
that the standard thin disc model is not self-consistent, and 
so attention has turned to other models such as advection- 
dominated flows (see Narayan 1997 for a review). 

Despite theoretical arguments that the thin disc is un- 
stable, thin disc spectra are widely used to fit observations 
of black hole candidates. Indeed, there is indirect evidence 
for (theoretically unstable) thin discs in that some observa- 
tions are well fit by thin disc spectra. For example, some 
galactic black hole candidates such as Nova Muscae have X- 



ray spectra that are well fit by multi-temperature thin disc 
models (R. Narayan, private communication). 

One can stabilize radiation pressure dominated discs 
by modifying the usual prescription for the shear stress 
t r( f, ~ ap so that t r ^, ~ ap g (p g = gas pressure) or some 
combination of gas and radiation pressure (e.g. Lightman & 
Eardley 1974, Piran 1978). Arguments have been advanced 
in favour of such a modification by Eardley & Lightman 
(1975), Coroniti (1981), Sakimoto & Coroniti (1981,1989), 
and Stella & Rosner (1984). These arguments rely on the 
thermodynamic peculiarities of magnetic buoyancy in a ra- 
diation dominated plasma. More recent work has vastly in- 
creased our understanding of magnetically driven angular 
momentum diffusion in discs (see the review of Balbus & 
Hawley 1997). Numerical experiments suggest that t r ^ is 
limited by Lorentz forces rather than by buoyant escape 
of magnetic fields (Stone et al. 1996), although these ex- 
periments do not include radiation pressure and radiative 
diffusion. 

An alternative route to viscous stability (e.g. Liang 
1977) is to modify the disc cooling law using convection, 
although this is not successful in eliminating the thermal in- 
stability (Piran 1978). Modification of the disc cooling rate 
has not seemed a promising approach. 

Recently, however, while investigating neutron star po- 
lar cap accretion, Klein & Arons (1989, 1991) noticed the 
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development of evacuated regions, or "photon bubbles"^ in 
their radiation hydrodynamics simulations. A subsequent 
linear analysis (Arons 1992) revealed that the source of the 
photon bubbles was an overstable mode present in radiation 
dominated, magnetized atmospheres. Earlier incarnations of 
the photon bubble instability in an unmagnetized plasma 
have turned out to be flawed; Marzec (1978) gives a full 
discussion of this point. Later numerical studies of the non- 
linear evolution of the photon bubble instability (Hsu et al. 
1997) showed that vertical transport of energy is enhanced 
in the nonlinear outcome. 

One is naturally led to inquire whether this same pho- 
ton bubble instability is present in the radiation dominated 
parts of accretion discs, and if so, what the consequences 
might be for disc structure. Unfortunately Arons's (1992) 
analysis is not immediately applicable to discs. It considers 
physical parameters relevant to neutron star polar cap accre- 
tion: a thermal timescale long compared to the dynamical 
timescale, negligible gas pressure, and superthermal mag- 
netic field. Our purpose in this paper is to generalize Arons's 
work to the regime appropriate to discs. 

First we evaluate conditions inside accretion discs us- 
ing a standard a model (§2). Then we write down a set of 
governing equations appropriate to these conditions (§3). In 
§4 we construct a model equilibrium to perturb about, and 
in §§5 and 6 we work out the linear theory in the WKB 
limit. The astrophysical implications of the result depend 
on the nonlinear outcome, about which we speculate in §7. 
§8 contains a summary. 



2 CONDITIONS IN ACCRETION DISCS 

What conditions are relevant to a study of the photon bub- 
ble instability in accretion discs? Consider accretion onto 
a black hole of mass m Mq at a rate rhL e ddC 2 , where 
L e dd = 4:TtGMc/K e3 , and K e s is the electron scattering opac- 
ity. Using the standard thin-disc, one-zone model (Shakura 
& Sunyaev 1973), assuming radiation pressure is dominant 
and k — k 



0.4 gem , we find that at radius rGM/c 



E ~ 0.37 a 1 m 1 r 3/ ' 2 gem 2 , 

T ~ 2.4 x 10 7 a- 1/4 m- 1/4 r- 3/8 K, 

a i in-7 -1 -1 • - 2 3/2 -3 

p ~ 4.2 x 10 a m m r ' gem 



(1) 
(2) 
(3) 



Here and throughout we ignore corrections due to general 
relativity and the inner boundary of the disc. Radiation- 
dominated accretion discs have entropy profiles that sug- 
gest convective instability (e.g. Bisnovatyi-Kogan & Blin- 
nikov 1977). Calculations of disc vertical structure using 
mixing-length theory (Shakura, Sunyaev, & Zilitinkevich 
1978) show, however, that the vertical radiative flux of en- 
ergy always dominates the vertical advective flux and so the 
mean structure of the disc is not likely to be very different 
from the values calculated above. 

The electron scattering optical depth is 



t s ~ Efi: es /2 = 0.074 a 



1 • -1 3/2 

m r ' 



(4) 



* Bubble is something of a misnomer, since surface tension plays 
no role. The phenomenon is really more like convection. 



while the "true" optical depth, which measures the extent 
to which radiation is thermalized inside the disc, depends 
on the absorption opacity K a : 



E^ K es K a jl. 



(5) 



The absorption opacity is a complicated function of density, 
temperature, metallicity and the radiation spectrum. Taking 
K a — Kp, the Planck mean opacity, and evaluating Kp from 
tables for a solar composition gas (Magee et al. 1995) we find 
that r* ranges from < 1 to ^> 1. For example, at r — 50 from 
a 10 Mq black hole accreting at m = 1, and taking a = 0.1, 
we find r* ~ 14.8. Discs around supermassive black holes 
have somewhat higher effective optical depths, since they 
are cooler and the Planck mean opacity increases sharply 
below 3 x 10 6 K due to bound-free absorption by metals. All 
this implies that the true optical depth of perturbations that 
are smaller than the scale height can be small, and so the 
radiation field must be treated using an approximation that 
is valid in this regime. 

The ratio of gas pressure p g to total pressure P is 



a — In i a w in-6 -1/4 -1/4 . -2 21/8 

p r = Pg/P = 1.6 x 10 a ' m ' m r ' 



(6) 



so the boundary of the radiation-dominated region lies at 

(7) 



-,nr\ 2/21 2/21 • 16/21 

r c ~ 160 a ' m m 1 



Also the sound speed c s is 
c s /c — 3.0 rhr~ 3//2 . 



(8) 



The radiative diffusivity is Do = c/(/tp). In dimensionless 
form, 

D 



M 



c a H 



4.5 a, 



(9) 



where H = disc scale height. Thus Mo ~ 1, while in neutron 
star polar cap accretion Mo <C 1. This difference changes 
the character of the instability significantly, and is the most 
important difference between our work and Arons's work. 
We also need an estimate for the magnetic field 



strength. We use the prescription (3m 



' /„,2 



(4a)" 1 

(cf = P/p is the isothermal sound speed; va = Alfven 
speed), consistent with the simulations of Hawley, Gammie, 
& Balbus (1995). So if a = 0.1, (3 M = 2.5. For reasonable 
values of a, then, /3m ~ 1. Arons (1992) focused on the 
case 13m <C 1, but developed a more general analysis in an 
appendix. 

Finally, heat conduction, radiative viscosity, ordinary 
viscosity, and ordinary resistivity are all completely negligi- 
ble. 

To summarize the results of this section, accretion 
discs have approximately thermal magnetic fields, so (3m = 
C i/ V A ~ 1. In their inner region radiation pressure domi- 
nates, so f3 r <g 1. The thermal timescale is comparable to the 
dynamical timescale (to within a factor of a) so the dimen- 
sionless radiative diffusion rate Mo = (c/ (up)) / (ciH) ~ 1. 
Finally, the thermalization length h = (p^/ftfta) -1 varies 
widely but can be ~ H. Since we will consider perturba- 
tions on scales small compared to the scale height, the per- 
turbations can also have a scale small compared to the ther- 
malization length. Thus it is not clear a priori that it is 
appropriate to use the Rosseland diffusion approximation 
for the radiation field. 



3 BASIC EQUATIONS 

Since we cannot use the Rosseland diffusion approximation 
for the radiation field, we turn to the more general flux- 
limited nonequilibrium diffusion approximation (see Mihalas 
& Mihalas 1984 and references therein) . This approximation 
is very similar to the Rosseland diffusion approximation (it 
still sets Kij = J8ij/3, where K and J are the second and 
zeroth angular moments of the intensity), but, loosely speak- 
ing, it allows for the possibility that the radiation field has 
a different "temperature" than the gas. More precisely, it 
does not require that J — oTg/iv, where T g = gas temper- 
ature. Because it is a diffusion approximation, it is strictly 
valid only on scales large compared to the photon mean free 
path (pre es ) -1 , but it gives qualitatively sensible results on 
smaller scales. 

The governing equations for the gas, then, are the con- 
tinuity equation, 

Ap = -p(V-v), (10) 
where D t = dt + (v • V), the momentum equation, 
(B ■ V)B VB 2 4ttkp 



pD t v = -Vp s - pV0 - 



4ir 



8. + ^ H ' < U > 



where H = J* Jn/(47r) ^| is the frequency-integrated flux 
and 4> is the gravitational potential, and the gas energy equa- 
tion (u = internal energy per unit volume), 



D t u = -7«(V • v) + 4irK a p(J - B), 



(12) 



where J is the frequency-integrated mean intensity, B ~ 
aTg/n, and we approximate the absorption opacity by the 
Planck mean opacity. The magnetic field evolution is gov- 
erned by the induction equation 



D t B = -B(V • v) + (B • V)v 



(13) 



and the constraint V ■ B = 0. The mean intensity evolution 
is given by 

-Dt J — ( Dtp = —V • H + K a p(B — J), (14) 
c \ Spc I 



and the flux evolution is given by 

-D t H = --VJ- npU, 
c 3 



(15) 



where we approximate k by the Rosseland mean opacity 



4 MODEL PROBLEM: STRATIFIED 
ATMOSPHERE 

Magnetized, radiation-dominated accretion discs are dy- 
namically evolving flows, since they are subject to the mag- 
netorotational instability (Balbus & Hawley 1991) and, pos- 
sibly, ordinary convective instability. Linear theory therefore 
cannot even in principle provide a rigorous guide to their dy- 
namics. The best we can hope for is to find a model problem 



that captures the essence of physical conditions in accretion 
discs and is simple enough to solve. 

The model problem we have chosen is a nonrotating 
stratified atmosphere with V(j> = gz = const., a uniform 
magnetic field, and a constant flux H z from below. This 
is the least complicated model that is potentially subject 
to the photon bubble instability. It allows us to focus on 
photon bubbles alone, disentangled from the magnetorota- 
tional instability, magnetic Rayleigh- Taylor instability, and 
convective instability. 

The model equilibrium is determined by the vertical 
momentum equation 



1 dp 4-KKp 
-tt H ti z 

p dz c 



9 



0, 



the gas energy equation 
J = B, 



(16) 



(17) 



the intensity equation 

V-H = 0, (18) 
and the vertical component of the flux equation 

ldJ_ 

3 dz 

This set of equations admits a two-parameter family of 
solutions. A natural set of parameters is j3 r evaluated at 
z — z , and C = H z 4n n / (eg) , which is the ratio of the flux 
to the critical flux where the effective gravity vanishes (i.e., 
the local Eddington limit). Defining the pressure scale height 



= —KpH z 



(19) 



H 



at z = z we have 



dz 

hp = - H ^ = 

az 

which imply 



C 



1 1 r 



-H 



-H 



d\nP r 

dz 

dlnP g 

dz 



4 1-/3/ 
± {-, £ (l-3fl./4) 

0r I (1 - Pr) 



1 - c 



(20) 
(21) 

(22) 
(23) 



The entropy profile suggests convective stability (en- 
tropy increases upwards) if 



C < C c 



4( 7 - l)(4-7/3 r + 3/3 r . 2 ) 



(24) 



" (7 - 1)(16 - 12/3 r - 3/3 r 2 ) + I3 r 2 

(Kutter 1970, Wentzel 1970, Tayler 1954). In the limit of 
small p r this amounts to 



C< l-/3 r + 0(/3 r 2 ). 



(25) 



For model atmospheres in which Mo 1 this implies con- 
vective stability. In discs the question of convective stability 
is more subtle, and outside the scope of this paper, because 
of the presence of radiative damping and rotation. 



t H denotes the scale height and Hi a component of the flux. 
Likewise B denotes the thermal mean intensity and Bi a com- 
ponent of the magnetic field. The difference should be clear in 
context. 



5 LINEAR THEORY 

We will consider only short-wavelength (WKB) perturba- 
tions. Longer-wavelength modes are global in nature and 
depend on the details of the equilibrium; but the careful 



construction of disc equilibria is a futile endeavour since 
all weakly magnetized discs are unstable. The perturbations 
have the form 



Sf ~ exp 



i(k x x + J k z {z')dz — uit) 



(26) 



where / is any one of the perturbed variables. Like Arons, 
we take k x = (1 — p 2 ) 1 ^ 2 k and k z = pk. The basic small 
parameter of the WKB approximation is {kH)~ x ~ e <JC 1 
(we assume that p ~ 1). 

Retaining all terms that are potentially of leading order, 
the linearized continuity equation is 



— lUJOp — OVz—rf- = —IpVL ■ OV 

H 

the momentum equation is 

—iupSv — —ik.Spg 

{SnpH + kS p H + np5H) - 

the induction equation is 

— icuSB = -iB(k • 5v) + i(k ■ B)<5v, 
the gas energy equation is 

— iljSu - &v z ^j- = — i7?i(k • <5v) + Ann a p{SJ — SB), (30) 
the intensity equation is 

li-iuSJ - 6v,!*2-) - |^Hpk • 5v) = 

-ik-5H + K a (SB-6J), (31) 

and the flux equation is 

- (iioSU) = —ikSJ - npSH- KSpH - SKpR. (32) 

c o 

Together with the constraint k ■ <5B = 0, these equa- 
tions imply a complicated eleventh-order dispersion relation 
D(u>, k) =0 which we shall not record here (the most com- 
pact form of the full relation is the above equations). We 
have confirmed that it contains the dispersion relations for 
magnetohydrodynamic (MHD) waves, magnetoatmospheric 
waves (e.g. Thomas 1982), internal waves, sound waves in a 
radiating fluid (Mihalas & Mihalas 1984), and the overstablc 
photon bubble mode of Arons (1992) as special cases. 



(27) 



gSpz, 



(28) 



(29) 



6 OVERSTABILITY 
6.1 Numerical Solution 

The full dispersion relation is analytically intractable. In 
the end we will need to expand it in a small parameter to 
retrieve the relevant pieces that describe the photon bubble 
mode. To motivate an asymptotic approach, we shall first 
solve the dispersion relation numerically. The parameters 
are those appropriate to a disc around a 10 Mq black hole 
accreting at the Eddington rate, at r = 50GM/c 2 . Assume 
a = 0.1, so /3m = 2.5, and that the field is purely vertical. 
Then j3 r = 0.05, c/a = 120, h = 0.07H, 7 = 5/3, M = 0.4. 
Take C — 1 — 3/3 r /2 = 0.93 so that convection is absent. We 
consider a set of modes with varying k and p — 1/3. The 
real part of the full dispersion relation is shown in Figure 1. 
The real part of the overstable mode is marked with a heavy 
solid line; the imaginary part has comparable magnitude. 



At large wavenumber the dispersion relation is easily 
interpreted because it becomes analytic. For general p, 

(c 2 -icV)(u, 2 -^c 2 fcV)x 

(c 4 - (/3 r7 c 2 fc 2 + -^c 2 k 2 )u 2 + J_/3 r7C 2 fcV) x (33) 

(a,+ g) a (a,+ 4M(1 -^r 1)Jo )=0. 



The first term in parentheses is the flux-limited electromag- 
netic wave with phase velocity c/y/3 (see Mihalas & Mihalas 
1984). The second term is the Alfven wave. The third term 
contains the fast and slow MHD modes. These are labeled in 
Figure 1. Notice that the overstable mode becomes the slow 
MHD mode at large wavenumber. The final three modes are 
strongly damped entropy modes. 

Notice that the slow MHD mode has rather low fre- 
quency. This is because at short wavelengths radiation dif- 
fuses rapidly out of the perturbation and so the radiation 
pressure perturbation is nil. Thus only gas pressure provides 
a restoring force for this mode. For /3 r C 1, the slow mode 
has u> 2 w (3 r c 2 fk 2 p 2 , so the slow mode velocity is directly 
related to the sound speed associated with the gas pressure 
alone. A radiation pressure dominated fluid is thus rather 
delicate on small scales in that it is easily compressed. 

6.2 Vertical Magnetic Field 

We have shown that an overstable mode exists; we will now 
demonstrate the existence of this mode analytically. Again 
we consider only the simple case of a purely vertical mag- 
netic field. The basic small parameter is e = (kH)^ 1 . The 
survey of conditions in accretion discs (§2) then suggests 
the following scalings for the other parameters in the prob- 
lem. We take H ~ l,c» ~ l,c;/c ~ e, /3m ~ 1, Z* ~ l,/3r ~ 

J3/2 



M 



1, C = 1 - 0{/3 r ), and p 



1. Using a little 

1/2 



asymptotic foresight, we also take u> ~ e 

Expanding D(uj,\s.) through leading order in e, we find 
the remarkably simple relation 



~igkp(l — p 2 



(34) 



One root of this equation, the one with negative phase ve- 
locity, describes the overstable photon bubble mode. 

What happens if we increase the importance of gas pres- 
sure? Suppose that (3 r ~ e. Then we find 



u! 2 = —igkp(l — p 2 ) + (3 r c 2 p 2 k 2 . 



(35) 



The first term is the photon bubble term, while the second 
is that of a slow MHD mode in which gas pressure provides 
the only restoring force. If we increase (3 r /{kH) still further, 
then the first term becomes subdominant. This branch of 
the dispersion relation then becomes the slow mode, which 
is stable to leading order in WKB. 

Since the overstability is no longer present to leading 
order in e when k > (/3 r iJ) _1 , and k > l/H, the oversta- 
bility fails in a WKB sense when /3 r ~ 1. This provides 
an approximate limit on the overstable region in parameter 
space. 

Now suppose that the diffusion time is long compared to 
a dynamical time, i.e. Mo ~ e 1 ^ 2 . Again expanding D(ui, k) 
to leading order in e we find 



2 , • 4Cj/i . 2 n 

uj + vjJ — ^ \-igkp(l — p ) 



0. 



(36) 



The new term causes damping. It shows that instability 
is present in a WKB sense only on scales such that k > 
l/(HMi). 

6.3 General Discussion 

We can obtain a better physical understanding of the origin 
of the photon bubble instability, and a more general dis- 
persion relation, by expanding the linearized equations in e. 
This is not trivial because we must assign an explicit relative 
ordering in e to all the perturbed variables. We do this by 
solving for the eigenvectors Sf(Sp) (/ represents any of the 
perturbed variables), and then determining the relative or- 
dering. We use the same ordering for the model parameters 
as that described at the beginning of §6.2. This procedure 
reveals that the photon bubble mode consists predominantly 
of motions along the magnetic field. This is easy to under- 
stand because the frequency of the mode, ~ \fgk~, is smaller 
by order e 1 / 2 than the Alfven frequency at the same scale, 
~ kvA- Thus the field is stiff enough to resist any motion 
that bends the field lines. 

Writing out the dominant terms in full and allowing for 
a general orientation of the initially uniform magnetic field, 
the continuity equation becomes 



iujSp = — ipk ■ Sv. 



(37) 



The perturbed velocities perpendicular to the magnetic field 
vanish. Denoting the unit vector parallel to the field by b, 
the parallel component of the velocity is governed by 



icjSv ■ b 



(38) 



so the material feels only the perturbed radiation force. The 
induction equation is irrelevant, since the field is stiff. The 
perturbed flux obeys 



k • 5H = 0. 



(39) 



The perturbed flux is conserved, so there is no exchange 
of energy between the radiation field and the fluid on a 
timescale Combining this result with the perturbed 

intensity equation, we find 



5H = -^ 



H k 



(k-H) 

k 2 



(40) 



In words, the fluid feels a radiation force that is directed 
parallel to the wave crests and is inversely proportional to 
the density perturbation. Radiation escapes more rapidly 
along density minima, while energy flow is impeded along 
density maxima. Combining eqns.(^), (pS8|), and (|40[), and 
denoting the unit vector along k by k, we find the general 
dispersion relation 

uj 2 = -ikg(k • b) [b • z - (k • z)(k ■ b)] . (41) 

This implies that the growth rate of the instability is largest 
for k z = (where the WKB approximation is not strictly 
valid) and for B at an angle of ±7r/4 to the wavevector. 

It is readily shown that if the fluid is allowed to move 
freely, unconstrained by the magnetic field, then the mode 
frequency vanishes if we start from the reduced set of equa- 
tions above. More rigorously, if we begin with the full set of 
linearized equations and turn off the magnetic field, then the 



stability criterion reduces to the condition that the Brunt- 
Vaisala frequency be real. 

To understand the overstability on a more qualitative 
level, consider the history of a single fluid element in the 
case k z = where B makes an angle 7r/4 with the vertical. 
Suppose the fluid element initially lies in a density mini- 
mum, so it feels an increased radiation force from below. It 
is accelerated upward along the field line. This soon puts it 
in a region of convergent flow and then in a density maxi- 
mum: the phase velocity of the overstable mode is such that 
the density maxima progress downward along field lines. It 
is now shadowed by the density maximum, the radiation 
flux drops somewhat, and it falls downward along the field 
line. This puts it in a region of divergent flow and then a 
density minimum, so it is accelerated upward again. Each 
round of acceleration is larger than the one before, and so 
an overstability results. 

The photon bubble instability described by Arons is 
physically similar to ours, but there are two differences. 
The most important is that Arons considers the asymptotic 
regime M <C 1 appropriate to neutron star polar caps. This 
leads to a rather more complicated dispersion relation for 
the overstable modes, with the growth rate proportional to 
Mq . In this limit, there is a small phase offset between the 
radiation flux and the density, rather than the large phase 
offset that initiates the instability in the limits M ~ 1. An- 
other difference is that Arons uses the Rosseland diffusion 
approximation, which requires that the gas and radiation 
have the same temperature, while we use a nonequilibrium 
diffusion approximation, which allows the gas and radiation 
"temperature" to differ. In the end this turns out to have 
no effect whatsoever on the overstable modes, because gas 
pressure (and hence gas temperature) is completely negli- 
gible for modes with frequency and wavelength comparable 
to the photon bubble mode. Finally, we note that in the 
appropriate limit we are able to recover Arons's dispersion 
relation, thus confirming his analysis. 



7 NONLINEAR OUTCOME 

The astrophysical implications of the photon bubble insta- 
bility depend on the nonlinear outcome, to which the linear 
theory is an unreliable guide. Numerical experiments by Hsu 
et al. (1997) show that in the regime Mo = c/(T es Ci) <C 1, 
Pm -C 1 the photon bubble instability leads to greatly en- 
hanced vertical transport of energy. It thus seems likely that 
the photon bubble instability will enhance vertical energy 
transport in discs. 

It is tempting to make an analogy between photon 
bubbles and convection. Convection generates such efficient 
transport of energy that it erases the inverted entropy gra- 
dient that initiated it and, if forced, maintains the connec- 
tive fluid in a marginally stable state. It is natural to think 
that the photon bubble instability drives the disc toward a 
marginally stable state as well. Our analysis does not give 
any rigorous stability criteria, but it does show that when 
p r ~ 1 the instability is no longer present at leading order. 
Thus the instability is likely to be absent or greatly reduced 
in strength when j3 r ~ 1. If the disc is initially radiation 
dominated, photon bubbles might then transport energy ef- 



ficiently out of the disc, lowering the temperature and raising 
f3 r toward 1. 

The nonlinear outcome of the photon bubble instability 
can be characterized by a cooling rate Q~(E, T c ), which is 
the energy lost per unit time per unit disc area. The ther- 
mal and viscous stability properties of the disc depend on 
how Q~ varies with temperature and surface density (Piran 
1978). Ultimately this can only be evaluated from a fully 
nonlinear theory, but if Q~ is a steep enough function of T c , 
then the disc can be viscously and thermally stable. 

Our analysis also shows that when the dimensionless 
radiative diffusion rate Mo decreases the instability weakens, 
in that instability is only present for kH > M ~ 2 . This effect 
may also shut off photon bubbles. At low accretion rates, 
for example, a disc with /3 r ~ 1 will have large Mo, while at 
larger accretion rates Mo is smaller. 

Photon bubbles might also change the emergent spec- 
trum. The bubbles make the disc porous to radiation, so a 
photon traverses a shorter path to the surface than it would 
if the disc were subject only to ordinary radiative difffu- 
sion. The photon distribution therefore has less opportunity 
to thermalize. In addition, the vertical transport might be 
more episodic, enhancing variability. 

This discussion is speculative. Other physical processes 
may contribute to vertical energy transport in discs: ordi- 
nary convection (but see the discussion of Rees 1987) or 
magnetic Rayleigh- Taylor instability (see the flux-tube cal- 
culation of Sakimoto & Coroniti 1989) may dominate pho- 
ton bubbles. In addition, MHD turbulence initiated by the 
Balbus-Hawley instability must coexist with photon bub- 
bles. We have developed a linear theory only when these 
effects are absent. The full nonlinear development of the 
radiation dominated disc can probably only be studied real- 
istically via three dimensional numerical experiments. Nu- 
merical methods exist for treating the radiative transfer in 
a flux-limited Rosseland diffusion approximation. Since the 
character of the unstable mode is identical in the nonequi- 
librium and Rosseland diffusion approximations, numerical 
studies of the nonlinear evolution of a radiation dominated 
disc may be immediately practical. 
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FIGURE CAPTIONS 

8 SUMMARY 

We have considered the linear theory of a radiation domi- 
nated atmosphere with spatially constant magnetic field as a 
model for the radiation dominated inner parts of thin accre- 
tion discs around compact objects. The model is subject to 
an overstable photon bubble mode that tends to separate ra- 
diation and matter. The photon bubble dispersion relation 
for a general orientation of the magnetic field is given by 
eqn.(Ulj). Vertical energy transport is likely to be enhanced 
in the nonlinear outcome. The disc may then deflate until it 
is no longer radiation pressure dominated. 
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Figure 1. Figure 1. The real part of the WKB dispersion relation 
for a radiation-dominated atmosphere with uniform vertical field. 
See text for details. The abscissa is the wavenumber in units of 
the scale height; the ordinate is the phase velocity in units of 
the isothermal sound speed (P/p) 1 / 2 . The real part of the phase 
velocity for the photon bubble mode is shown as a heavy line. 
The growth rate (imaginary part) is comparable in magnitude. 



